From 4e017862e159f75e43a0d254564b5e9c2c6c0b37 Mon Sep 17 00:00:00 2001 From: tsteven4 <13596209+tsteven4@users.noreply.github.com> Date: Tue, 13 Aug 2024 14:21:31 -0600 Subject: [PATCH] use implicit conversions for lat,lon degs<->rads in grtcirc interface (#1321) * use implict conversions for lat,lon degs<->rads. * use Position* return values. * convert grtcirc i/f to PositionDeg, PositionRad. * fix mac compile, denote converting ctor. * inline PositionDeg/Rad converting ctors. * respect PosotionX encapsulation. * differentiate member names for PositionDeg, PositiionRad. This can catch unintended usages where the type is not what was expected. * clean up unneccesary usage of Waypoint::position. * delete obsolete declaration. * spell * update garmin include * consistenly unpack Positions to consts in grtcirc. * provide default initializers for PositionX. --- arcdist.cc | 27 +++------- arcdist.h | 3 +- bend.cc | 26 +++------ defs.h | 45 ++++++++++++++++ garmin.cc | 1 - gdb.cc | 4 +- grtcirc.cc | 140 ++++++++++++++++++++++++------------------------- grtcirc.h | 37 +++++-------- igc.cc | 3 +- interpolate.cc | 17 ++---- kml.cc | 3 +- position.cc | 6 +-- position.h | 4 -- radius.cc | 4 +- radius.h | 5 -- route.cc | 8 +-- smplrout.cc | 33 +++++------- trackfilter.cc | 43 +++++---------- waypt.cc | 5 +- xcsv.cc | 14 ++--- xcsv.h | 3 +- 21 files changed, 194 insertions(+), 237 deletions(-) diff --git a/arcdist.cc b/arcdist.cc index 5f56b8ca4..d213222d6 100644 --- a/arcdist.cc +++ b/arcdist.cc @@ -44,8 +44,7 @@ void ArcDistanceFilter::arcdist_arc_disp_wpt_cb(const Waypoint* arcpt2) { static const Waypoint* arcpt1 = nullptr; - double prjlat; - double prjlon; + PositionDeg prjpos; double frac; if (arcpt2 && arcpt2->latitude != BADVAL && arcpt2->longitude != BADVAL && @@ -62,12 +61,8 @@ void ArcDistanceFilter::arcdist_arc_disp_wpt_cb(const Waypoint* arcpt2) if (ed->distance == BADVAL || projectopt || ed->distance >= pos_dist) { double dist; if (ptsopt) { - dist = gcdist(RAD(arcpt2->latitude), - RAD(arcpt2->longitude), - RAD(waypointp->latitude), - RAD(waypointp->longitude)); - prjlat = arcpt2->latitude; - prjlon = arcpt2->longitude; + dist = gcdist(arcpt2->position(), waypointp->position()); + prjpos = arcpt2->position(); frac = 1.0; } else { if (waypointp == nullptr) { @@ -77,13 +72,9 @@ void ArcDistanceFilter::arcdist_arc_disp_wpt_cb(const Waypoint* arcpt2) fatal(FatalMsg() << "Internal error: Attempt to project waypoint without predecessor"); } - dist = linedistprj(arcpt1->latitude, - arcpt1->longitude, - arcpt2->latitude, - arcpt2->longitude, - waypointp->latitude, - waypointp->longitude, - &prjlat, &prjlon, &frac); + std::tie(dist, prjpos, frac) = linedistprj(arcpt1->position(), + arcpt2->position(), + waypointp->position()); } /* convert radians to float point statute miles */ @@ -92,8 +83,7 @@ void ArcDistanceFilter::arcdist_arc_disp_wpt_cb(const Waypoint* arcpt2) if (ed->distance > dist) { ed->distance = dist; if (projectopt) { - ed->prjlatitude = prjlat; - ed->prjlongitude = prjlon; + ed->prjpos = prjpos; ed->frac = frac; ed->arcpt1 = arcpt1; ed->arcpt2 = arcpt2; @@ -171,8 +161,7 @@ void ArcDistanceFilter::process() wp->wpt_flags.marked_for_deletion = 1; removed++; } else if (projectopt) { - wp->longitude = ed->prjlongitude; - wp->latitude = ed->prjlatitude; + wp->SetPosition(ed->prjpos); if (!arcfileopt && (ed->arcpt2->altitude != unknown_alt) && (ptsopt || (ed->arcpt1->altitude != unknown_alt))) { diff --git a/arcdist.h b/arcdist.h index 2458ffc9e..7e30e4e70 100644 --- a/arcdist.h +++ b/arcdist.h @@ -44,8 +44,7 @@ private: struct extra_data { double distance; - double prjlatitude; - double prjlongitude; + PositionDeg prjpos; double frac; const Waypoint* arcpt1; const Waypoint* arcpt2; diff --git a/bend.cc b/bend.cc index 774bea92e..09325673f 100644 --- a/bend.cc +++ b/bend.cc @@ -55,8 +55,8 @@ void BendFilter::init() Waypoint* BendFilter::create_wpt_dest(const Waypoint* wpt_orig, const Waypoint* wpt_orig_adj) const { - double distance = radtometers(gcdist(RAD(wpt_orig->latitude), RAD(wpt_orig->longitude), - RAD(wpt_orig_adj->latitude), RAD(wpt_orig_adj->longitude))); + double distance = radtometers(gcdist(wpt_orig->position(), + wpt_orig_adj->position())); if (distance <= maxDist) { return nullptr; } @@ -64,10 +64,7 @@ Waypoint* BendFilter::create_wpt_dest(const Waypoint* wpt_orig, const Waypoint* double frac = maxDist / distance; auto* wpt_dest = new Waypoint(*wpt_orig); - linepart(wpt_orig->latitude, wpt_orig->longitude, - wpt_orig_adj->latitude, wpt_orig_adj->longitude, - frac, - &wpt_dest->latitude, &wpt_dest->longitude); + wpt_dest->SetPosition(linepart(wpt_orig->position(), wpt_orig_adj->position(), frac)); return wpt_dest; } @@ -76,19 +73,10 @@ int BendFilter::is_small_angle(const Waypoint* wpt_orig, const Waypoint* wpt_orig_prev, const Waypoint* wpt_orig_next) const { - double lat_orig = RAD(wpt_orig->latitude); - double long_orig = RAD(wpt_orig->longitude); - - double lat_orig_prev = RAD(wpt_orig_prev->latitude); - double long_orig_prev = RAD(wpt_orig_prev->longitude); - - double lat_orig_next = RAD(wpt_orig_next->latitude); - double long_orig_next = RAD(wpt_orig_next->longitude); - - double heading_prev = heading_true_degrees(lat_orig, long_orig, - lat_orig_prev, long_orig_prev); - double heading_next = heading_true_degrees(lat_orig, long_orig, - lat_orig_next, long_orig_next); + double heading_prev = heading_true_degrees(wpt_orig->position(), + wpt_orig_prev->position()); + double heading_next = heading_true_degrees(wpt_orig->position(), + wpt_orig_next->position()); double heading_diff = heading_next - heading_prev; diff --git a/defs.h b/defs.h index df3e8c3aa..906bd9a61 100644 --- a/defs.h +++ b/defs.h @@ -19,10 +19,12 @@ #ifndef DEFS_H_INCLUDED_ #define DEFS_H_INCLUDED_ +#include // for nan #include // for NULL, nullptr_t, size_t #include // for int32_t, uint32_t #include // for NULL, fprintf, FILE, stdout #include // for time_t +#include // for inv_pi, pi #include // for optional #include // for move @@ -107,6 +109,14 @@ constexpr double MPH_TO_MPS(double a) { return a * kMPSPerMPH;} /* knots(nautical miles/hour) to meters/second */ constexpr double KNOTS_TO_MPS(double a) {return a * kMPSPerKnot;} +/* Degrees to radians */ +constexpr double kDegreesPerRadian = 180.0 * std::numbers::inv_pi; +constexpr double DEG(double x) { return x * kDegreesPerRadian; } + +/* Radians to degrees */ +constexpr double kRadiansPerDegree = 1.0 / kDegreesPerRadian; +constexpr double RAD(double x) { return x * kRadiansPerDegree; } + constexpr int kDatumOSGB36 = 86; // GPS_Lookup_Datum_Index("OSGB36") constexpr int kDatumWGS84 = 118; // GPS_Lookup_Datum_Index("WGS 84") @@ -246,6 +256,35 @@ struct bounds { double min_alt; /* -unknown_alt => invalid */ }; +struct PositionRad; // forward declare +struct PositionDeg +{ + PositionDeg() = default; + explicit(false) inline PositionDeg(const PositionRad& posr); /* converting ctor */ + PositionDeg(double latd, double lond) : latD(latd), lonD(lond) {} + + double latD{std::nan("")}; + double lonD{std::nan("")}; +}; + +struct PositionRad +{ + PositionRad() = default; + explicit(false) inline PositionRad(const PositionDeg& posd); /* converting ctor */ + PositionRad(double latr, double lonr) : latR(latr), lonR(lonr) {} + + double latR{std::nan("")}; + double lonR{std::nan("")}; +}; + +inline PositionDeg::PositionDeg(const PositionRad& posr) : + latD(posr.latR * kDegreesPerRadian), + lonD(posr.lonR * kDegreesPerRadian) {} + +inline PositionRad::PositionRad(const PositionDeg& posd) : + latR(posd.latD * kRadiansPerDegree), + lonR(posd.lonD * kRadiansPerDegree) {} + /* * This is a waypoint, as stored in the GPSR. It tries to not * cater to any specific model or protocol. Anything that needs to @@ -318,6 +357,12 @@ public: void SetCreationTime(qint64 t, qint64 ms = 0); Geocache* AllocGCData(); int EmptyGCData() const; + PositionDeg position() const {return PositionDeg(latitude, longitude);} + void SetPosition(const PositionDeg& pos) + { + latitude = pos.latD; + longitude = pos.lonD; + } // mimic std::optional interface, but use our more space // efficient wp_flags. diff --git a/garmin.cc b/garmin.cc index 7c20fc836..13012f5ee 100644 --- a/garmin.cc +++ b/garmin.cc @@ -43,7 +43,6 @@ #include "garmin_fs.h" // for garmin_fs_garmin_after_read, garmin_fs_garmin_before_write #include "garmin_tables.h" // for gt_find_icon_number_from_desc, PCX, gt_find_desc_from_icon_number #include "geocache.h" // for Geocache, Geocache::type_t, Geocache... -#include "grtcirc.h" // for DEG #include "jeeps/gpsapp.h" // for GPS_Set_Baud_Rate, GPS_Init, GPS_Pre... #include "jeeps/gpscom.h" // for GPS_Command_Get_Lap, GPS_Command_Get... #include "jeeps/gpsmem.h" // for GPS_Track_Del, GPS_Way_Del, GPS_Pvt_Del diff --git a/gdb.cc b/gdb.cc index 47f0d072a..ca276ce86 100644 --- a/gdb.cc +++ b/gdb.cc @@ -258,9 +258,7 @@ GdbFormat::gdb_add_route_waypt(route_head* rte, Waypoint* ref, const int wpt_cla /* At this point we have found a waypoint with same name, but probably from another data stream. Check coordinates! */ - double dist = radtometers(gcdist( - RAD(ref->latitude), RAD(ref->longitude), - RAD(tmp->latitude), RAD(tmp->longitude))); + double dist = radtometers(gcdist(ref->position(), tmp->position())); if (fabs(dist) > 100) { fatal(MYNAME ": Route point mismatch!\n" \ diff --git a/grtcirc.cc b/grtcirc.cc index 90c0e0850..b8e4e34cb 100644 --- a/grtcirc.cc +++ b/grtcirc.cc @@ -25,9 +25,10 @@ #include // for errno, EDOM #include // for cos, sin, fabs, atan2, sqrt, asin, atan, isnan #include // for pi -#include // for make_tuple, tuple +#include // for tie, tuple, make_tuple, ignore -#include "defs.h" // for METERS_TO_MILES +#include "defs.h" // for PositionRad, DEG, METERS_TO_MILES, PositionDeg +#include "grtcirc.h" static constexpr double EARTH_RAD = 6378137.0; @@ -46,17 +47,6 @@ static double dotproduct(double x1, double y1, double z1, return (x1 * x2 + y1 * y2 + z1 * z2); } -/* - * Note: this conversion to miles uses the WGS84 value for the radius of - * the earth at the equator. - * (radius in meters)*(100cm/m) -> (radius in cm) - * (radius in cm) / (2.54 cm/in) -> (radius in in) - * (radius in in) / (12 in/ft) -> (radius in ft) - * (radius in ft) / (5280 ft/mi) -> (radius in mi) - * If the compiler is half-decent, it'll do all the math for us at compile - * time, so why not leave the expression human-readable? - */ - double radtomiles(double rads) { constexpr double radmiles = METERS_TO_MILES(EARTH_RAD); @@ -68,10 +58,15 @@ double radtometers(double rads) return (rads * EARTH_RAD); } -double gcdist(double lat1, double lon1, double lat2, double lon2) +double gcdist(PositionRad pos1, PositionRad pos2) { errno = 0; + const double lat1 = pos1.latR; + const double lon1 = pos1.lonR; + const double lat2 = pos2.latR; + const double lon2 = pos2.lonR; + double sdlat = sin((lat1 - lat2) / 2.0); double sdlon = sin((lon1 - lon2) / 2.0); @@ -92,8 +87,13 @@ double gcdist(double lat1, double lon1, double lat2, double lon2) /* This value is the heading you'd leave point 1 at to arrive at point 2. * Inputs and outputs are in radians. */ -static double heading(double lat1, double lon1, double lat2, double lon2) +static double heading(PositionRad pos1, PositionRad pos2) { + const double lat1 = pos1.latR; + const double lon1 = pos1.lonR; + const double lat2 = pos2.latR; + const double lon2 = pos2.lonR; + double v1 = sin(lon2 - lon1) * cos(lat2); double v2 = cos(lat1) * sin(lat2) - sin(lat1) * cos(lat2) * cos(lon2 - lon1); /* rounding error protection */ @@ -106,10 +106,9 @@ static double heading(double lat1, double lon1, double lat2, double lon2) return atan2(v1, v2); } -/* As above, but outputs is in degrees from 0 - 359. Inputs are still radians. */ -double heading_true_degrees(double lat1, double lon1, double lat2, double lon2) +double heading_true_degrees(PositionRad pos1, PositionRad pos2) { - double h = 360.0 + DEG(heading(lat1, lon1, lat2, lon2)); + double h = 360.0 + DEG(heading(pos1, pos2)); if (h >= 360.0) { h -= 360.0; } @@ -120,11 +119,9 @@ double heading_true_degrees(double lat1, double lon1, double lat2, double lon2) // Note: This is probably not going to vectorize as it uses statics internally, // so it's hard for the optimizer to prove it's a pure function with no side // effects, right? -double linedistprj(double lat1, double lon1, - double lat2, double lon2, - double lat3, double lon3, - double* prjlat, double* prjlon, - double* frac) +std::tuple linedistprj(PositionRad pos1, + PositionRad pos2, + PositionRad pos3) { static double _lat1 = -9999; static double _lat2 = -9999; @@ -144,17 +141,22 @@ double linedistprj(double lat1, double lon1, double dot; - *prjlat = lat1; - *prjlon = lon1; - *frac = 0; + /* prjpos must be of type PositionRad as we assign + * prjpos.lat, prjpos.lon below in radians. + */ + PositionRad prjpos = pos1; - /* degrees to radians */ - lat1 = RAD(lat1); - lon1 = RAD(lon1); - lat2 = RAD(lat2); - lon2 = RAD(lon2); - lat3 = RAD(lat3); - lon3 = RAD(lon3); + double frac = 0; + + /* we use these values below assuming they are in radians, + * => posn must be of type PositionRad. + */ + const double lat1 = pos1.latR; + const double lon1 = pos1.lonR; + const double lat2 = pos2.latR; + const double lon2 = pos2.lonR; + const double lat3 = pos3.latR; + const double lon3 = pos3.lonR; bool newpoints = true; if (lat1 == _lat1 && lat2 == _lat2 && lon1 == _lon1 && lon2 == _lon2) { @@ -183,10 +185,7 @@ double linedistprj(double lat1, double lon1, /* 'a' is the axis; the line that passes through the center of the earth * and is perpendicular to the great circle through point 1 and point 2 * It is computed by taking the cross product of the '1' and '2' vectors.*/ - auto [xt, yt, zt] = crossproduct(x1, y1, z1, x2, y2, z2); - xa = xt; - ya = yt; - za = zt; + std::tie(xa, ya, za) = crossproduct(x1, y1, z1, x2, y2, z2); la = sqrt(xa * xa + ya * ya + za * za); if (la) { @@ -230,15 +229,15 @@ double linedistprj(double lat1, double lon1, * atan2 would be overkill because lp and fabs(dot) are both * known to be positive. */ - *prjlat = DEG(asin(yp)); + prjpos.latR = asin(yp); if (xp == 0 && zp == 0) { - *prjlon = 0; + prjpos.lonR = 0; } else { - *prjlon = DEG(atan2(zp, xp)); + prjpos.lonR = atan2(zp, xp); } - *frac = d1 / (d1 + d2); + frac = d1 / (d1 + d2); - return atan(fabs(dot) / lp); + return {atan(fabs(dot) / lp), prjpos , frac}; } /* otherwise, get the distance from the closest endpoint */ @@ -264,34 +263,32 @@ double linedistprj(double lat1, double lon1, } if (fabs(d1) < fabs(d2)) { - return gcdist(lat1, lon1, lat3, lon3); + return {gcdist(pos1, pos3), prjpos, frac}; } else { - *prjlat = DEG(lat2); - *prjlon = DEG(lon2); - *frac = 1.0; - return gcdist(lat2, lon2, lat3, lon3); + prjpos = pos2; + frac = 1.0; + return {gcdist(pos2, pos3), prjpos, frac}; } } else { /* lp is 0 when 3 is 90 degrees from the great circle */ - return std::numbers::pi / 2; + return {std::numbers::pi / 2, prjpos, frac}; } } else { /* la is 0 when 1 and 2 are either the same point or 180 degrees apart */ dot = dotproduct(x1, y1, z1, x2, y2, z2); if (dot >= 0) { - return gcdist(lat1, lon1, lat3, lon3); + return {gcdist(pos1, pos3), prjpos, frac}; } else { - return 0; + return {0, prjpos, frac}; } } } -double linedist(double lat1, double lon1, - double lat2, double lon2, - double lat3, double lon3) +double linedist(PositionRad pos1, PositionRad pos2, PositionRad pos3) { - double dummy; - return linedistprj(lat1, lon1, lat2, lon2, lat3, lon3, &dummy, &dummy, &dummy); + double dist; + std::tie(dist, std::ignore, std::ignore) = linedistprj(pos1, pos2, pos3); + return dist; } /* @@ -301,20 +298,20 @@ double linedist(double lat1, double lon1, * Ref: http://mathworld.wolfram.com/RotationFormula.html */ -void linepart(double lat1, double lon1, - double lat2, double lon2, - double frac, - double* reslat, double* reslon) +PositionDeg linepart(PositionRad pos1, PositionRad pos2, double frac) { - /* result must be in degrees */ - *reslat = lat1; - *reslon = lon1; - - /* degrees to radians */ - lat1 = RAD(lat1); - lon1 = RAD(lon1); - lat2 = RAD(lat2); - lon2 = RAD(lon2); + /* respos must be of type PositionRad as we assign + * respos.lat, respos.lon below in radians. + */ + PositionRad respos = pos1; + + /* we use these values below assuming they are in radians, + * => posn must be of type PositionRad. + */ + const double lat1 = pos1.latR; + const double lon1 = pos1.lonR; + const double lat2 = pos2.latR; + const double lon2 = pos2.lonR; /* polar to ECEF rectangular */ double x1 = cos(lon1) * cos(lat1); @@ -357,11 +354,12 @@ void linepart(double lat1, double lon1, yr = std::clamp(yr, -1.0, 1.0); zr = std::clamp(zr, -1.0, 1.0); - *reslat = DEG(asin(yr)); + respos.latR = asin(yr); if (xr == 0 && zr == 0) { - *reslon = 0; + respos.lonR = 0; } else { - *reslon = DEG(atan2(zr, xr)); + respos.lonR = atan2(zr, xr); } } + return respos; } diff --git a/grtcirc.h b/grtcirc.h index 981f87228..9d8900dd0 100644 --- a/grtcirc.h +++ b/grtcirc.h @@ -22,35 +22,24 @@ #ifndef GRTCIRC_H #define GRTCIRC_H -#include // for inv_pi +#include // for tuple +#include "defs.h" // for PositionRad, PositionDeg -double gcdist(double lat1, double lon1, double lat2, double lon2); -double heading_true_degrees(double lat1, double lon1, double lat2, double lon2); +/* Note PositionDeg and PositionRad can be implicity converted to + * each other, so you may use either to interface to these functions. + */ + +double gcdist(PositionRad pos1, PositionRad pos2); +double heading_true_degrees(PositionRad pos1, PositionRad pos2); -double linedistprj(double lat1, double lon1, - double lat2, double lon2, - double lat3, double lon3, - double* prjlat, double* prjlon, - double* frac); +std::tuple linedistprj(PositionRad pos1, + PositionRad pos2, + PositionRad pos3); -double linedist(double lat1, double lon1, - double lat2, double lon2, - double lat3, double lon3); +double linedist(PositionRad pos1, PositionRad pos2, PositionRad pos3); double radtometers(double rads); double radtomiles(double rads); -void linepart(double lat1, double lon1, - double lat2, double lon2, - double frac, - double* reslat, double* reslon); - -/* Degrees to radians */ -constexpr double kDegreesPerRadian = 180.0 * std::numbers::inv_pi; -constexpr double DEG(double x) { return x * kDegreesPerRadian; } - -/* Radians to degrees */ -constexpr double kRadiansPerDegree = 1.0 / kDegreesPerRadian; -constexpr double RAD(double x) { return x * kRadiansPerDegree; } - +PositionDeg linepart(PositionRad pos1, PositionRad pos2, double frac); #endif diff --git a/igc.cc b/igc.cc index 23b026470..f93bb793c 100644 --- a/igc.cc +++ b/igc.cc @@ -869,8 +869,7 @@ int IgcFormat::correlate_tracks(const route_head* pres_track, const route_head* // Get a crude indication of groundspeed from the change in lat/lon int deltat_msec = (*wpt_rit)->GetCreationTime().msecsTo(wpt->GetCreationTime()); speed = (deltat_msec == 0) ? 0: - radtometers(gcdist(RAD(wpt->latitude), RAD(wpt->longitude), - RAD((*wpt_rit)->latitude), RAD((*wpt_rit)->longitude))) / + radtometers(gcdist(wpt->position(), (*wpt_rit)->position())) / (0.001 * deltat_msec); if (global_opts.debug_level >= 2) { printf(MYNAME ": speed=%.2fm/s\n", speed); diff --git a/interpolate.cc b/interpolate.cc index f1367c1f2..4a88ac0ba 100644 --- a/interpolate.cc +++ b/interpolate.cc @@ -65,8 +65,7 @@ void InterpolateFilter::process_rte(route_head* rte) } // And add them back, with interpolated points interspersed. - double lat1 = 0; - double lon1 = 0; + PositionDeg pos1; double altitude1 = unknown_alt; gpsbabel::DateTime time1; bool first = true; @@ -92,10 +91,7 @@ void InterpolateFilter::process_rte(route_head* rte) // interpolate even if time is running backwards. npts = std::abs(*timespan) / max_time_step; } else if (opt_dist != nullptr) { - double distspan = radtomiles(gcdist(RAD(lat1), - RAD(lon1), - RAD(wpt->latitude), - RAD(wpt->longitude))); + double distspan = radtomiles(gcdist(pos1, wpt->position())); npts = distspan / max_dist_step; } if (!std::isfinite(npts) || (npts >= INT_MAX)) { @@ -119,11 +115,7 @@ void InterpolateFilter::process_rte(route_head* rte) } else { wpt_new->creation_time = gpsbabel::DateTime(); } - linepart(lat1, lon1, - wpt->latitude, wpt->longitude, - frac, - &wpt_new->latitude, - &wpt_new->longitude); + wpt_new->SetPosition(linepart(pos1, wpt->position(), frac)); if (altspan.has_value()) { wpt_new->altitude = altitude1 + (frac * *altspan); } else { @@ -142,8 +134,7 @@ void InterpolateFilter::process_rte(route_head* rte) track_add_wpt(rte, wpt); } - lat1 = wpt->latitude; - lon1 = wpt->longitude; + pos1 = wpt->position(); altitude1 = wpt->altitude; time1 = wpt->creation_time.toUTC(); // use utc to avoid tz conversions. } diff --git a/kml.cc b/kml.cc index cacf0385c..ae27864b0 100644 --- a/kml.cc +++ b/kml.cc @@ -2070,8 +2070,7 @@ void KmlFormat::wr_position(Waypoint* wpt) } else { Waypoint* newest_posn= posn_trk_head->waypoint_list.back(); - if (radtometers(gcdist(RAD(wpt->latitude), RAD(wpt->longitude), - RAD(newest_posn->latitude), RAD(newest_posn->longitude))) > 50) { + if (radtometers(gcdist(wpt->position(), newest_posn->position())) > 50) { track_add_wpt(posn_trk_head, new Waypoint(*wpt)); } else { /* If we haven't move more than our threshold, pretend diff --git a/position.cc b/position.cc index a7275b84e..759e9e11f 100644 --- a/position.cc +++ b/position.cc @@ -49,10 +49,8 @@ void PositionFilter::position_runqueue(const WaypointList& waypt_list, int qtype for (int j = i + 1 ; j < nelems ; ++j) { if (!qlist.at(j).deleted) { - double dist = gc_distance(qlist.at(j).wpt->latitude, - qlist.at(j).wpt->longitude, - qlist.at(i).wpt->latitude, - qlist.at(i).wpt->longitude); + double dist = radtometers(gcdist(qlist.at(j).wpt->position(), + qlist.at(i).wpt->position())); if (dist <= pos_dist) { if (check_time) { diff --git a/position.h b/position.h index 8c1f8111f..4c2204172 100644 --- a/position.h +++ b/position.h @@ -57,10 +57,6 @@ private: /* Member Functions */ - static double gc_distance(double lat1, double lon1, double lat2, double lon2) - { - return radtometers(gcdist(RAD(lat1), RAD(lon1), RAD(lat2), RAD(lon2))); - } void position_runqueue(const WaypointList& waypt_list, int qtype); /* Data Members */ diff --git a/radius.cc b/radius.cc index 5552790d2..7e978970b 100644 --- a/radius.cc +++ b/radius.cc @@ -35,8 +35,8 @@ void RadiusFilter::process() { foreach (Waypoint* waypointp, *global_waypoint_list) { - double dist = gc_distance(waypointp->latitude, waypointp->longitude, - home_pos->latitude, home_pos->longitude); + double dist = radtomiles(gcdist(waypointp->position(), + home_pos->position())); if ((dist >= pos_dist) == (exclopt == nullptr)) { waypointp->wpt_flags.marked_for_deletion = 1; diff --git a/radius.h b/radius.h index 8bebd63f9..3e18953c9 100644 --- a/radius.h +++ b/radius.h @@ -51,11 +51,6 @@ private: /* Member Functions */ - static double gc_distance(double lat1, double lon1, double lat2, double lon2) - { - return radtomiles(gcdist(RAD(lat1), RAD(lon1), RAD(lat2), RAD(lon2))); - } - /* Data Members */ double pos_dist{}; diff --git a/route.cc b/route.cc index c4c25f56e..34b728a47 100644 --- a/route.cc +++ b/route.cc @@ -278,16 +278,12 @@ computed_trkdata track_recompute(const route_head* trk) /* * gcdist and heading want radians, not degrees. */ - double tlat = RAD(thisw->latitude); - double tlon = RAD(thisw->longitude); - double plat = RAD(prev->latitude); - double plon = RAD(prev->longitude); if (!thisw->course_has_value()) { // Only recompute course if the waypoint // didn't already have a course. - thisw->set_course(heading_true_degrees(plat, plon, tlat, tlon)); + thisw->set_course(heading_true_degrees(prev->position(), thisw->position())); } - double dist = radtometers(gcdist(plat, plon, tlat, tlon)); + double dist = radtometers(gcdist(prev->position(), thisw->position())); tdata.distance_meters += dist; /* diff --git a/smplrout.cc b/smplrout.cc index 6a71295b3..eb5388126 100644 --- a/smplrout.cc +++ b/smplrout.cc @@ -97,18 +97,15 @@ double SimplifyRouteFilter::compute_track_error(const neighborhood& nb) const switch (metric) { case metric_t::crosstrack: track_error = radtomiles(linedist( - wpt1->latitude, wpt1->longitude, - wpt2->latitude, wpt2->longitude, - wpt3->latitude, wpt3->longitude)); + wpt1->position(), + wpt2->position(), + wpt3->position())); break; case metric_t::length: track_error = radtomiles( - gcdist(RAD(wpt1->latitude), RAD(wpt1->longitude), - RAD(wpt3->latitude), RAD(wpt3->longitude)) + - gcdist(RAD(wpt3->latitude), RAD(wpt3->longitude), - RAD(wpt2->latitude), RAD(wpt2->longitude)) - - gcdist(RAD(wpt1->latitude), RAD(wpt1->longitude), - RAD(wpt2->latitude), RAD(wpt2->longitude))); + gcdist(wpt1->position(), wpt3->position()) + + gcdist(wpt3->position(), wpt2->position()) - + gcdist(wpt1->position(), wpt2->position())); break; case metric_t::relative: default: // eliminate false positive warning with g++ 11.3.0: ‘error’ may be used uninitialized in this function [-Wmaybe-uninitialized] @@ -119,19 +116,15 @@ double SimplifyRouteFilter::compute_track_error(const neighborhood& nb) const (wpt1->GetCreationTime() != wpt2->GetCreationTime())) { double frac = static_cast(wpt1->GetCreationTime().msecsTo(wpt3->GetCreationTime())) / static_cast(wpt1->GetCreationTime().msecsTo(wpt2->GetCreationTime())); - double reslat; - double reslon; - linepart(wpt1->latitude, wpt1->longitude, - wpt2->latitude, wpt2->longitude, - frac, &reslat, &reslon); - track_error = radtometers(gcdist( - RAD(wpt3->latitude), RAD(wpt3->longitude), - RAD(reslat), RAD(reslon))); + auto respos = linepart(wpt1->position(), + wpt2->position(), + frac); + track_error = radtometers(gcdist(wpt3->position(), respos)); } else { // else distance to connecting line track_error = radtometers(linedist( - wpt1->latitude, wpt1->longitude, - wpt2->latitude, wpt2->longitude, - wpt3->latitude, wpt3->longitude)); + wpt1->position(), + wpt2->position(), + wpt3->position())); } // error relative to horizontal precision track_error /= (6 * wpt3->hdop); diff --git a/trackfilter.cc b/trackfilter.cc index bcb7b77d0..0d43a42b4 100644 --- a/trackfilter.cc +++ b/trackfilter.cc @@ -531,8 +531,7 @@ void TrackFilter::trackfilter_split() if (distance > 0) { double curdist = radtometers( - gcdist(RAD(prev_wpt->latitude), RAD(prev_wpt->longitude), - RAD(wpt->latitude), RAD(wpt->longitude))); + gcdist(prev_wpt->position(), wpt->position())); if (curdist <= distance) { new_track_flag = false; } else if constexpr(TRACKF_DBG) { @@ -600,10 +599,8 @@ void TrackFilter::trackfilter_move() void TrackFilter::trackfilter_synth() { - double last_course_lat; - double last_course_lon; - double last_speed_lat = std::nan(""); /* Quiet gcc 7.3.0 -Wmaybe-uninitialized */ - double last_speed_lon = std::nan(""); /* Quiet gcc 7.3.0 -Wmaybe-uninitialized */ + PositionDeg last_course_pos; + PositionDeg last_speed_pos; gpsbabel::DateTime last_speed_time; int nsats = 0; @@ -626,18 +623,14 @@ void TrackFilter::trackfilter_synth() wpt->reset_speed(); } first = false; - last_course_lat = wpt->latitude; - last_course_lon = wpt->longitude; - last_speed_lat = wpt->latitude; - last_speed_lon = wpt->longitude; + last_course_pos = wpt->position(); + last_speed_pos = wpt->position(); last_speed_time = wpt->GetCreationTime(); } else { if (opt_course) { - wpt->set_course(heading_true_degrees(RAD(last_course_lat), - RAD(last_course_lon),RAD(wpt->latitude), - RAD(wpt->longitude))); - last_course_lat = wpt->latitude; - last_course_lon = wpt->longitude; + wpt->set_course(heading_true_degrees(last_course_pos, + wpt->position())); + last_course_pos = wpt->position(); } if (opt_speed) { if (last_speed_time.msecsTo(wpt->GetCreationTime()) != 0) { @@ -649,14 +642,11 @@ void TrackFilter::trackfilter_synth() // Note that points with the same time can occur because the input // has truncated times, or because we are truncating times with // toTime_t(). - wpt->set_speed(radtometers(gcdist( - RAD(last_speed_lat), RAD(last_speed_lon), - RAD(wpt->latitude), - RAD(wpt->longitude))) / + wpt->set_speed(radtometers(gcdist(last_speed_pos, wpt->position())) + / (0.001 * std::abs(last_speed_time.msecsTo(wpt->GetCreationTime()))) ); - last_speed_lat = wpt->latitude; - last_speed_lon = wpt->longitude; + last_speed_pos = wpt->position(); last_speed_time = wpt->GetCreationTime(); } else { wpt->reset_speed(); @@ -888,10 +878,7 @@ void TrackFilter::trackfilter_faketime() bool TrackFilter::trackfilter_points_are_same(const Waypoint* wpta, const Waypoint* wptb) { return - radtometers(gcdist(RAD(wpta->latitude), - RAD(wpta->longitude), - RAD(wptb->latitude), - RAD(wptb->longitude))) < kDistanceLimit && + radtometers(gcdist(wpta->position(), wptb->position())) < kDistanceLimit && std::abs(wpta->altitude - wptb->altitude) < 20 && wpta->courses_equal(*wptb) && wpta->speeds_equal(*wptb) && @@ -909,10 +896,8 @@ void TrackFilter::trackfilter_segment_head(const route_head* rte) for (auto it = wptlist.cbegin(); it != wptlist.cend(); ++it) { auto* wpt = *it; if (it != wptlist.cbegin()) { - double cur_dist = radtometers(gcdist(RAD(prev_wpt->latitude), - RAD(prev_wpt->longitude), - RAD(wpt->latitude), - RAD(wpt->longitude))); + double cur_dist = radtometers(gcdist(prev_wpt->position(), + wpt->position())); // Denoise points that are on top of each other, // keeping the first and last of the group. diff --git a/waypt.cc b/waypt.cc index e7c8c1898..b26455511 100644 --- a/waypt.cc +++ b/waypt.cc @@ -200,11 +200,12 @@ waypt_add_url(Waypoint* wpt, const QString& link, const QString& url_link_text, wpt->AddUrlLink(UrlLink(link, url_link_text, url_link_type)); } +// TODO: change inputs to PositionDeg? double gcgeodist(const double lat1, const double lon1, const double lat2, const double lon2) { - return radtometers(gcdist(RAD(lat1), RAD(lon1), RAD(lat2), RAD(lon2))); + return radtometers(gcdist(PositionDeg(lat1, lon1), PositionDeg(lat2, lon2))); } /* @@ -355,7 +356,7 @@ double waypt_course(const Waypoint* A, const Waypoint* B) { if (A && B) { - return heading_true_degrees(RAD(A->latitude), RAD(A->longitude), RAD(B->latitude), RAD(B->longitude)); + return heading_true_degrees(A->position(), B->position()); } else { return 0; } diff --git a/xcsv.cc b/xcsv.cc index c64d9c4d8..afc3cf373 100644 --- a/xcsv.cc +++ b/xcsv.cc @@ -963,8 +963,7 @@ void XcsvFormat::xcsv_resetpathlen(const route_head* head) { pathdist = 0; - oldlat = 999; - oldlon = 999; + old_position.reset(); csv_route = csv_track = nullptr; switch (xcsv_style->datatype) { case trkdata: @@ -993,12 +992,13 @@ XcsvFormat::xcsv_waypt_pr(const Waypoint* wpt) double utmn; char utmzc; - if (oldlon < 900) { - pathdist += radtometers(gcdist(RAD(oldlat),RAD(oldlon), - RAD(wpt->latitude),RAD(wpt->longitude))); + if (old_position) { + pathdist += radtometers(gcdist(old_position.value(), + wpt->position())); } - longitude = oldlon = wpt->longitude; - latitude = oldlat = wpt->latitude; + old_position = wpt->position(); + latitude = wpt->latitude; + longitude = wpt->longitude; QString write_delimiter; if (xcsv_style->field_delimiter == u"\\w") { diff --git a/xcsv.h b/xcsv.h index 838c60c07..b7d4467f8 100644 --- a/xcsv.h +++ b/xcsv.h @@ -366,8 +366,7 @@ private: XcsvFile* xcsv_file{nullptr}; const XcsvStyle* xcsv_style{nullptr}; double pathdist = 0; - double oldlon = 999; - double oldlat = 999; + std::optional old_position; int waypt_out_count = 0; const route_head* csv_track = nullptr; -- 2.30.2